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Abstract 

This paper presents a numerical compression strategy for the boundary integral equation of 
acoustic scattering in two dimensions. These equations have oscillatory kernels that we represent 
in a basis of wave atoms, and compress by thresholding the small coefficients to zero. 
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Left: a kite-shaped scatterer. Middle: a depiction of the kernel of the double- layer potential 
for this scatterer, sampled as a 1024x1024 matrix. Right: a zoomed-in view of the sparsity 
pattern of the "nonstandard wave atom matrix" , representing this kernel accurately using only 
60,000 matrix elements. 

This phenomenon was perhaps first observed in 1993 by Bradie, Coifman, and Grossman, 
in the context of local Fourier bases [5] . Their results have since then been extended in various 
ways. The purpose of this paper is to bridge a theoretical gap and prove that a well-chosen 
fixed expansion, the nonstandard wave atom form, provides a compression of the acoustic single 
and double layer potentials with wave number k as 0(k)-by-0(k) matrices with 0(fc 1+1 /°°) 
nonnegligible entries, with a constant that depends on the relative £2 accuracy e in an acceptable 
way. The argument assumes smooth, separated, and not necessarily convex scatterers in two 
dimensions. The essential features of wave atoms that enable to write this result as a theorem 
is a sharp time-frequency localization that wavelet packets do not obey, and a parabolic scaling 
wavelength ~ (essential diameter) 2 . Numerical experiments support the estimate and show that 
this wave atom representation may be of interest for applications where the same scattering 
problem needs to be solved for many boundary conditions, for example, the computation of 
radar cross sections. 



Acknowledgements. The first author is partially supported by an NSF grant. The second 
author is partially supported by an NSF grant, a Sloan Research Fellowship, and a startup grant 
from the University of Texas at Austin. 



1 



1 Introduction 



This paper is concerned with the sparse representation of the boundary integral operator of two- 
dimensional scattering problems. Let D be a bounded soft scatterer in M 2 with a smooth boundary 
and u mc {x) be the incoming wave field. The scattered field u{x) satisfies the two-dimensional 
exterior Dirichlet problem of the Helmholtz equation: 

-Au(x) - k 2 u{x) = in R d \ D, 

u{x) = -u mc {x) for x £ dD, 

^ lim \x\ 1/2 (Jj^ > Vn(x)^j - iku(x)j = 0. 

One attractive method for dealing with this problem is to reformulate it using a boundary integral 
equation for an unknown field (f)(x) on dD: 



\m + - ivG(x,y)^j <P(y)dy = -u mc {x), (1) 

where n y stands for the exterior normal direction of dD at the point y, and rj is a coupling constant of 
>(* 

equation and its normal derivative, given by 



order 0{k). The kernels G{x, y) and 9G q^' v ^ are, respectively, the Green's function of the Helmholtz 



G(x 1 y)= i -H^\k\\x-y\\ 



and 

■ )n -( x ,y) = ^HY > {k\\x-y 



dG , , ik (i) 



dn y y ' w/ 4 1 y " ||x-y|| y 

Once </>(x) is obtained from solving the integral equation, the scattered field u(x) at x £ M. 2 \D can 
be evaluated as 

u ( x ) = [ ( 9G ( X,y ^ - ir/G(x, y) J <^(y)dy. 



V dn y 

An important property of (1) from the computational point of view is that its condition number 
is often quite small and, as a result, one can advantageously solve (1) with an iterative algorithm 
like GMRES. At each step of the iterative solver, we need to apply the integral operator to a given 
function. Since the integral operator is dense, applying the operator directly is too expensive. In 
this paper, we address this issue by efficiently representing the operator as a sparse matrix in a 
system of wave atoms. 

Local cosines or wavelet packets have already been proposed for this task with great practical 
success, see Section 1.4 for some references, but we believe that the following two reasons make a 
case for wave atom frames: 

• The proposed construction is non-adaptive: wave atom frames of 1? are not designed for a 
specific value of k, and no optimization algorithm is needed to find a provably good basis. To 
achieve this result, the essential property of wave atoms is a parabolic scaling that we discuss 
later. 

• The choice of numerical realization for wave atoms follows some of the experience garnered 
throughout the 1990s in the study of local Fourier bases and wavelet packets. In particular, 
wave atoms offer a clean multiscale structure in the sense that they avoid the "frequency 
leakage" associated with wavelet packets defined from filterbanks. These aspects are discussed 
in [11]. 
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Of course, any non-adaptive all-purpose numerical compression method is likely to lag in per- 
formance behind an adaptive strategy that would include at least the former in its scope; but this 
is no excuse for discarding their study. Proper insight about architectures and scalings is important 
for designing the solution around which, for instance, a library of bases should be deployed for a 
best basis search. 

The main result of this paper says that the wave atom frame is in some sense near-optimal for 
representing the integral operator in (1) as a sparse matrix. Namely, full matrices would involve 
0(k 2 ) elements but we show that O(k) matrix elements suffice to represent G, and 0(k 1+6 ) matrix 
elements suffice to represent dG/dn y to a given accuracy e, for arbitrarily small 5 > 0. We believe 
that these bounds would not hold for wavelet packets, for instance, even if the best decomposition 
tree is chosen. In particular, wavelets would obviously not be suited for the job. 

Once the sparse representation of (1) is constructed in the wave atom frame, applying the 
operator requires only two wave atom transforms and has complexity 0(N 1+S ) for a problem with 
N unknowns (TV is proportional to k, 5 is arbitrarily small.) The implicit constant is in practice 
very small, and the wave atom transform itself has complexity 0(N log N). 

Such a strategy is attractive when the scattering problem needs to be solved several times 
with different incoming waves u mc {x). One important example is the computation of bistatic 
cross sections, where one needs to calculate the far field patterns of scattered fields for all possible 
incoming plane waves. 

1.1 Wave atoms 

Frames of wave atoms were introduced in [11] on the basis that they provide sparse representations 
of certain oscillatory patterns. As alluded to earlier, they are a special kind of oriented wavelet 
packets that do not suffer from the frequency leaking associated to filterbanks, and which obey the 
important parabolic scaling relation 

wavelength ~ (essential diameter) 2 . 

Let us rehearse the construction of wave atoms, and refer the reader to [11] for more details. In 
one dimension, wave atoms are an orthonormal basis indexed by the triple of integers A = (j, m, n). 
The construction is in the frequency domain; our convention for the Fourier transform is 

m= f e-^f(x)dx, f( x ) = -J—[ e™f(u)du. 

Jr" {^r Jr™ 

• First, j > is a scale parameter that should be thought of as indexing dilations of a factor 4 ', 
in other words, one should consider a first partition of the positive frequency axis into intervals 
of the form [c\2 23 , C22 23+2 \ (for some constants ci,C2 that will accommodate overlapping of 
basis functions). Choosing j such that frequency to is proportional to 2 2j is in contrast with 
wavelet theory, where to ~ 2 3 over the positive frequency support of a wavelet. 

• The parameter m with c\2 3 < m < C22 3+2 then indexes the further partitioning of each 
interval [ci2 2jl , C22 2j,+2 ] into 0(2 3 ) subintervals of size 0(2 3 ). More precisely, wave atoms are 
centered in frequency near ±uj\ where 

lu\ fa ir2 j m, c\2 3 < m < c 2 2 j+2 

and are compactly supported in the union of two intervals of length 2ir x 2 3 . The parabolic 
scaling is now apparent; the size of the support in frequency (~ 2 3 ) is proportional to the 
square root of the offset from the origin (~ 2 2j ). 
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• The parameter n £ Z indexes translations. A wave atom is centered in space near 

x\ = 2~ 3 n, 

and has essential support as narrow as the uncertainty principle allows, i.e., of length 0(2~- J ). 
We define to be the set of all admissible indices, i.e., 

n = {(j, m , n) : j > 0, a2 j <m< c 2 2 j+2 , n € Z}. 
Basis functions are then written 

<p x (x) = V/ 2 <p (j>m) @x-n)e™ mx , AGfi (2) 

where <P(j,m) depends weakly on j and m, and needs to be chosen adequately to form an orthobasis. 
The underlying delicate construction of the (fu m \ is due to Lars Villemoes [25] and summarized in 
[11]. 

In two dimensions, wave atoms are individually, but not collectively, formed as tensor products 
of the one-dimensional basis functions. The construction is "multiresolution" in the sense that 
there is only one dilation parameter; the indexing 5-uple of integers is fi = (j, m, n) where m = 
(mi, 777.2) and n = (721,722). More precisely, at scale j, the valid values for m = (mi, 777,2) satisfy 
< mi, 777,2 < C22 J+2 and Ci2 j < max(mi, 777,2). 

Wave atoms come as an orthonormal basis in two dimensions, but can be made fully directional — 
supported in a narrow cone in frequency with apex at the origin [2] — at the expense of increasing 
the redundancy to 2 or 4. The definition of such variants makes use of a unitary recombination 
involving Hilbert-transformed basis functions as in the definition of complex wavelet transforms, 
and is fully explained in [11]. 

None of the results of this paper would depend on the choice of variant; and for convenience we 
use the frame of wave atoms with redundancy four. With x = (x\,X2), the only property of wave 
atoms that we will need is the characterization 

^(x) = 2V( j ,m)(2^ 1 - ni, 2^2 - n 2 )e iVm - x , (3) 

where f(j : m) is a C°° non-oscillatory bump that depends on j and m, but in a non-essential manner, 
i.e., 

|<9> - m) (x)| < C aM (l + ||x||)- M , VM > 0, (4) 

with C at M independent of j and m. In addition, each tp(j t m) is simply the tensor product of two 
corresponding bumps for the ID transform. 

Although they may not necessarily form an orthonormal basis, wave atoms still form a tight 
frame in the sense that expanding a function is an isometry from L 2 (M?) to £2(1^), 

11/111 = Ekwi 2 

which is equivalent to 

f = ^2{f,<P»)<Pn- (5) 

The same properties hold in one dimension. 

The closest analogue to a "continuous wave atom transform" was introduced in the mathemati- 
cal literature by Cordoba and Fefferman in [8] . Wave atoms can be compared to brushlets [22] , but 
Villemoes's construction uses "local complex exponentials" instead of local cosines in frequency. 
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Discretized wave atoms are described in [11, 10]; they inherit the localization and tight-frame 
properties of their continuous counterpart. In particular, their bandlimited character confers an 
immediate control over the accuracy of computing inner products via quadrature. They come with 
fast FFT-based 0(N log N) algorithms for both the forward and adjoint transforms (see [11] for 
details). 

1.2 Operator expansions 

As functions can be analyzed and synthesized using coefficients, operators can also be expanded 
from matrix elements in a tight frame. When the frame is chosen appropriately, the expansion 
of the operator is approximately sparse: only a small percentage of the expansion coefficients are 
nonnegligible while all the rest of the coefficients are almost zero. Thresholding the small coefficients 
below a certain accuracy level gives rise to a sparse representation of the operator. Typically, this 
can be done in two ways. 

• The standard form of an operator A in the wave atom frame <p\ is 

xenx'en 

where ^4a,a' = (Atp\>, fx). For each fixed A, we define S(X) to be the set of all A' such that 
the modulus of A\^\/ is above a certain threshold. The sparse representation of A then takes 
the form 

x x'es(x) 

In practice, the sums in A and A' are also truncated in scale to account for the finite number 
of samples N of the functions to which A is applied. Since the wave atom transform computes 
the coefficients in the tight frame ip\ in 0(N log N) steps, the above equation naturally gives 
rise to an efficient method of applying the operator A to a given function /: 

— Apply the forward wave atom transform to compute the coefficients f\> := (f,<P\')- 

— For each A, compute g\ := Ylx'eS(X) A\ : yfy. The number of operations of this step is 
dictated by the number of nonnegligible coefficients of the representation of A. 

— Apply the adjoint wave atom transform to g\ to synthesize Af, i.e., Af Ylx fxdx- 

• The nonstandard form of A in the two dimensional frame (p^ is the set of coefficients A^ = 
J" R2 A(xi, X2)<p fJi (xi, X2) dx±dx2, such that the distributional kernel of A is expanded as 

A(xi,x 2 ) = ^ A M ^ M (xi,x 2 ) 

For a fixed threshold value, we define S to be the set of all fj, such that A^ is above the threshold 
in modulus. The sparse representation of A is now A(x\, X2) ~ J2fxeS AuW( x i> X2 )- Applying 
A to a given function / efficiently using this expansion is more involved than the case of the 
standard form. For a fixed index \i = ( j, m, n) with m = (771-1,7712) and n = (711,772), we 
define the two ID wave atom indices: 

Ai = (j, 777.1, Tii) and A^ = (j, m 2 , n 2 ). (6) 
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Since the two dimensional index m = (mi, 7712) satisfies < 7711,771,2 < c 2 2^ +2 and c\P < 
max(mi,m2), the set of all possible choices for A^ 1 and A2 are 

< mi < c 2 2 j+2 , and < m 2 < c 2 2 j+2 . (7) 

Some of the indices in (7) are not admissible, i.e., they are not part of the set of indices for the 
ID wave atom transform, since if fx = (j, m, n) corresponds to a ID wave atom it would need 
to satisfy ci2 J < m < c 2 2^ +2 . Non-admissible indices correspond to Gabor-type wave forms 
that partition the frequency domain uniformly and, hence, violate the parabolic scaling. We 
use £l e to denote this extended index set 

n e = {(j, m, n) : j > 0, < m < c 2 2 j+2 , n G Z}. 

The frame formed by these bases functions ip\ with A G VL e are called the extended wave atom 
frame. For a given function /, the computation of all the coefficients (/, (p\) with A G £l e can 
be done easily by extending the existing forward wave atom transform to include the extra 
771 indices (0 < m < ci2 J ) in 17 e . This forward extended wave atom transform still has an 
0(N log N) complexity since the number of extra non-admissible indices m at each scale j is 
a small fixed fraction of the number of existing m in £1 at scale j. The adjoint transform can 
be extended similarly and the resulting transform is called the adjoint extended wave atom 
transform. We note that these new transforms are not orthonormal any more since parts of 
the input functions are analyzed redundantly. 

Using the notation in (6) and the tensor-product property of ip^, we have 

A(x u x 2 ) ~ ^v4 m ^ a m(xi)^ a m(x 2 ). (8) 

When A is applied to a given function /, we have 

Af{xi) w J^^ a m(xi) f / (p x »(x 2 )f(x 2 )dx 2 j = (p x ^ ( / tp x ^(x 2 )f(x 2 )dx 2 

Using the extended transforms, we can derive from the above equation a fast algorithm for 
applying A to f using the nonstandard form: 

- Apply the forward extended wave atom transform to compute the coefficients fy := 
(f,ipy) for all indices A' G SI 6 . 

- For each A G S7 6 , compute g\ := X^gss.t. a m =a^m/a^- The number of operations in 
this step is proportional to the number of nonnegligible coefficients in the nonstandard 
expansion of the operator A. 

- Apply the inverse extended wave atom transform to synthesize Af from g\, i.e., Af ~ 

We would like to point out that nonstandard expansions only exist for two-dimensional frames 
that have a tensor product representation for each basis function, since the decomposition 
in (8) is essential for the derivation. For example, there is no known nonstandard form 
representation for an operator in the tight frame of curvelets [6] . 
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In what follows we focus exclusively on the nonstandard form, because of its relative simplicity 
over the standard form. This claim may seem paradoxical in view of the preceding discussion, 
but the basis functions <^(xi,X2) of the nonstandard form have a single dilation parameter and 
hence a fixed isotropic aspect ratio, unlike the tensor basis (px(xi)ip\>(x2) which come in all shapes 
and can be very elongated when j ^ f . The isotropy of the envelope of tp^ makes some of the 
stationary-phase argument in the sequel, simpler in our view. The kernel A(x±, X2) does not display 
the anisotropy of tensor wave atoms along the axes — if anything, A(x±,X2) is singular along the 
diagonal X\ = X2 — so choosing the nonstandard form also seems to us more natural from the 
numerical viewpoint. We can however not exclude at this point that the standard form may enjoy 
comparable sparsity properties as the nonstandard form. 

1.3 Sparsity of the nonstandard wave atom matrix 

In this section, we formulate the main result on sparsity of the nonstandard wave atom matrix of 
the acoustic single and double-layer potentials, in two dimensions. 

The scatterer is a union of closed, nonintersecting C°° curves = Ua=i ^ a em bedded in R 2 . 
For each a, assume that x(t) : I a 1— ► £l a is a C°° periodic parametrization of £l a , and take I a = [0, 1] 
for simplicity. 

We assume the following mild geometric regularity condition on the scatterer: there exists D > 
such that 

||x(s) -x(t)|| > D\e 27ris -e 27Tit \, (9) 

essentially meaning that the curve £l a defining the scatterer cannot intersect itself. We write 
d(s,t) = \e 2nts — e 2mt \ for the Euclidean distance on the unit circle. 
When s € I a , t G Ip with 1 < cc, (3 < n, let 

G ( S ,t) = ^ 1) (fc||x( S )-x( i )||)||x(t)||, (10) 

and 

Gi(a,t) = jH^ikMs) - x(t)||) jj^^ljj • "x(0 (11) 
The nonstandard wave atom matrices of Go and G\, restricted to a couple of intervals I a x Ip, 

are 

K° = (G ,^), Kl = {G u ^). 

Our main result below concerns the existence of e-approximants K® and K^, corresponding to the 
restriction of \i to sets Ao and Ai, i.e., with a = 0, 1, 

M \ otherwise, 

and chosen by definition such that 

\\K a -K a \\ hM <e, a = 0,1. (12) 

The £2 norm of a nonstandard wave atom matrix is equivalent to a Hilbert-Schmidt norm for the 
corresponding operator, by the tight-frame property of wave atoms. This norm is of course stronger 
than the operator L 2 -to-L 2 norm; and much stronger than the norm used in [5]. 

In what follows the notation A < B means A < C B for some constant C that depends only 
on the nonessential parameters. Similarly, the notation A < e _1 /°° means A < C M e~ 1/M for all 
M > 0. The constants may change from line to line. 
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Theorem 1. Assume the scatterer is smooth and geometrically regular in the sense of (9). In 
the notations just introduced, there exist sets Ao and Ai that define e-approximants of K® and K\ 
respectively, and whose cardinality obeys 



|A | < C" M 



(i) 



-j\ 2+1/M 



|Ai| < C\ 



M 



k l+l/M e -l/M + ^ 



\ 2/3+1/M 



(13) 



(14) 



/or all M > 0, and where G M , Cm depend only on M and the geometry of the scatterers. 



The terms proportional to k or fc 1+1 /°° are due to the oscillations when s / t, and the terms 
e -2-1 / 00 and (A;/e) 2 / 3+1 /°° are due to the kernels' singularities on the diagonal s = t. While the 
growth rate of k for the oscillations term is smaller for Go than for Gi, the growth rate of the 
diagonal contribution is smaller for G\ than for Go since 

fc 2/3 e -2/3 = fc l/3 fe l/3 e -2/3 < fc> £ - 2) < 3^ + -2y 

Theorem 1 can also be formulated using relative errors instead of absolute errors. This viewpoint 
is important for considering the composed kernel 

G( ,i)(s, t) = Gi(s, t) - irjG (s, t), 77 x k. 

Call K^' 1 ^ the nonstandard wave atom matrix of G(q ; i). In the following result, we quantify the 
number of terms needed to obtain the relative error estimate 

\\K a - K a \\ i2M < e\\K a \\ hM , a = 0,1, or (0, 1). 

Corollary 2. Let rj x k. In the assumptions and notations of Theorem 1, let Ao, Ai and A( 0) i) 
be the sets of wave atom coefficients needed to represent the operators Go, G±, resp. G( 0) i) up to 
relative accuracy e. Then 

(15) 

■ A .l+l/Af e -l/M + (A . e -2 ) l/3+l/Mj j (16) 

(17) 

for all M > 0, and where the constants depend only on M and on the geometry of the scatterers. 

We do not know if factors such as k 1 / 00 and e _1// °° could be replaced by log factors. 

The proofs of Theorem 1 and Corollary 2 occupy Section 2. The main ingredients are sparsity 
estimates in £ p , stationary phase considerations, vaguelette-type estimates adapted to wave atoms, 
and £2 correspondence scale-by-scale with wavelets. In Section 3 we present some numerical exper- 
iments that support the theory and establish wave atoms as a practical tool for solving scattering 
problems. 



I An 



<G° 



ill 



(ke 



-2\1+1/M 



Ail < G M 



A 



(0,1) 



< G 



(0,1) 

M 



(ke 



-2\1+1/M 
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1.4 Related work 



There has been a lot of work on sparsifying the integral operator of (1), or some variants of it, in 
appropriate bases. In [5], Bradie et al. showed that the operator becomes sparse in a local cosine 
basis. They proved that the number of coefficients with absolute value greater than any fixed e is 
bounded by O(klogk) when the constant depends on e. Notice that our result in Theorem 1 is 
stronger as the £2 norm is used instead in (12). In [3], Averbuch et al. extended the work in [5] by 
performing best basis search in the class of adaptive hierarchical local cosine bases. 

Besides the local cosine transform, adaptive wavelet packets have been used to sparsify the inte- 
gral operator as well. Deng and Ling [13] applied the best basis algorithm to the integral operator 
to choose the right one dimensional wavelet packet basis. Golik [17] independently proposed to ap- 
ply the best basis algorithm on the right hand side of the integral equation (1). Shortly afterwards, 
Deng and Ling [14] gave similar results by using a predefined wavelet packet basis that refines the 
frequency domain near k. All of these approaches work with the standard form expansion of the 
integral operator. Recently in [18], Huybrechs and Vandewalle used the best basis algorithm for two 
dimensional wavelet packets to construct a nonstandard sparse expansion of the integral operator. 
In all of these results, the numbers of nonnegligible coefficients in the expansions were reported to 
scale like <9(/c 4 / 3 ). However, our result shows that, by using the nonstandard form based on wave 
atoms, the number of significant coefficients scales like 0(/c 1+1 /°°). 

Most of the approaches on sparsifying (1) in well-chosen bases require the construction of the 
full integral operator. Since this step itself takes 0(k 2 ) operations, it poses computational difficulty 
for large k values. In [4], Beylkin et al. proposed a solution to the related problem of sparsifying the 
boundary integral operator of the Laplace equation. They successfully avoided the construction of 
the full integral operator by predicting the location of the large coefficients and applying a special 
one-point quadrature rule to compute the coefficients. The corresponding solution for the integral 
operator of the Helmholtz equation is still missing. 

There has been a different class of methods, initiated by Rokhlin in [23, 24], that requires no 
construction of the integral operator and takes 0(A;logA;) operations in 2D to apply the integral 
operator. A common feature of these methods [7, 15, 16, 23, 24] is that they partition the spatial 
domain hierarchically with a tree structure and compute the interaction between the tree nodes in 
a multiscale fashion: Whenever two nodes of the tree are well-separated, the interaction (of the 
integral operator) between them is either accelerated either by Fourier transform-type techniques 
[7, 23, 24] or by directional low rank representations [15, 16]. 

A criticism of the methods in [7, 15, 16, 23, 24] is that the constant in front of the complexity 
O(klogk) is often quite high. On the other hand, since the FFT-based wave atom transforms are 
extremely efficient, applying the operator in the wave atom frame has a very small constant once 
the nonstandard sparse representation is constructed. Therefore, for applications where one needs 
to solve the same Helmholtz equation with many different right hand sides, the current approach 
based on the wave atom basis can offer an competitive alternative. As mentioned earlier, one 
important example is the computation of the radar cross section. 

2 Sparsity analysis 

This section contains the proof of Theorem 1. The overarching strategy is to reduce the £2 approx- 
imation problem to an estimate of £ p sparsity through a basic result of approximation theory, the 
direct "Jackson" estimate 

H^-i^lh <C|A|H \\KX, 
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where H-f^Hp = Here refers to the approximation of where only the |A| largest 

terms in magnitude are kept, and the others put to zero. The inequality is valid for all values of 
< p < 2 for which ||K M || P is finite. For a proof, see [21], p. 390. 

If the £2 error is to be made less than e, it is enough to have in some i v space, < p < 2, 
and take the number of terms defining to be 

|A| >C p e^ WKpllf* (18) 

for some adequate C p > 0. The sequence will be split into several fragments that will be studied 
independently. For each of these fragments F in fi space, the inequality (18) will be complemented 
by an estimate of the form (f) < C p k q(j} \ for all p > po- Three scenarios will occur in the 

sequel: 

• If 

p = and q(p) = - - -, (19) 

then \F\ < fce _1// °°, which is the first term in (13). 

• If 

P0 = 1 and q(p)=0, (20) 
then |F| < e -2-i/oo^ w ] 1 j c ] 1 j s ^ e second term in (13). 

• If 

po = - and q(p) = --1 + 5, (21) 
z p 

for arbitrarily small 8 > 0, then \F\ < (fc/e) 2 / 3+1 /°°, which is the second term in (14). 

The problem is therefore reduced to identifying contributions in the sequence that obey 
either of the three estimates above. In what follows we focus on the kernel K = G$. We mention 
in Section 2.9 how the proof needs to be modified to treat the kernel G\. 

2.1 Smoothness of Hankel functions 

Bessel and Hankel functions have well-known asymptotic expansions near the origin and near 
infinity. That these asymptotic behaviors also determine smoothness in a sharp way over the whole 
half-line is perhaps less well-known, so we formulate these results as lemmas that we prove in the 
Appendix. 

Lemma 1. For every integers m > and n > there exists C m ^ n > such that for all k > 0, 

f d \m ( C mn {kx)- l ' 2 x- m ifkx>l; 

I f — J \e- lkx H£\kx)\ \ <l C mn (kx)- n x~ m if0<kx<landm + n> 0; (22) 
^ ' \ C (1 + I log kx\) if < kx < 1 and m = n = 0. 

The same results hold if 1 is replaced by any number c > in kx < 1 vs. kx > 1. 

The point of equation (22) is that C min is independent of k. Slightly more regularity can be 
obtained near the origin when multiplying with the adequate power of x, as the following lemma 
shows in the case 
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Lemma 2. For every integer m > there exists C m > such that, for < x < 1, 

{C m if m = 0,1; 

C 2 (l + |logx|) i/m = 2; (23) 
C m x 2 ~ m ifm>2. 

Finally, we will need the following lower bound. 

Lemma 3. For each n > 0, i/iere exist c n > and C n > such that, when x > c n , 

\H^(x)\ > CnX- 1 ' 2 . 

2.2 Dyadic partitioning 

Consider K as in Theorem 1, with s € I a and t e Ip. If a = (3, K presents a singularity on its 
diagonal, whereas if a ^ (5 it presents no such singularity. The case a = /3 is representative and is 
treated in the sequel without loss of generality. 

In this section we assume, as we have above, that I a = [0, 1]. The first step of the proof is to 
partition the periodized square I a x I a at each scale j, into dyadic squares denoted 

Q = [2-^i,2-% + 1)] x [2-^,2^(92 + !)]■ 

We define wq a window localized near Q through 

wq(s, t) = w^s - qi,2H - q 2 ), 

where w is compactly supported on [—1, 2] 2 and of class C°°. As a result wq is compactly supported 
inside 

3Q = [2-% - l),2-''( gi + 2)] x [2^(92 " 1), ^(q 2 + 2)]. 

We also write xq = (2~iqi, 2~^q 2 ) for the bottom-left corner of Q, not to be confused with x, which 
is in physical space. 

Denote by Qj the set of dyadic squares at scale j; we assume that w is chosen so that we have 
the scale-by-scale partition of unity property 

w q = 1 - 

The kernel is now analyzed at each scale j as K = ^2q £ q. Kq, where Kq = wqK. Dyadic 
squares can be classified according to their location with respect to the diagonal s = t, where K is 
singular. 

1. Diagonal squares. Dyadic squares will be considered "diagonal squares" as soon as the dis- 
tance from their center to the diagonal s = t is less than 1/k. Scale-by-scale, this condition 
reads d{2~^q\, 2~^q 2 ) < 3max(2~ J , |). (We need to use of the circle distance d since q\ and 
q 2 are defined modulo 2K) There are 0{2 3 : max(l, 2^/k)) such diagonal squares at scale j. 
They correspond to the case kx < 1 in Lemma 1. 

2. Nondiagonal squares. When d(2~^qi, 2~iq 2 ) > 3max(2^, 1/k), we say the square is non- 
diagonal, in those squares, the kernel Kq is C°° but oscillatory. There are 0(2 2 - J ) such 
nondiagonal squares at scale j. They correspond to the case kx > 1 in Lemma 1. 



n 



We take the scale j of the dyadic partitioning to match the scale j in the wave atom expansion; 
namely if fi = (j, m, n), then 

When < p < 1, an estimate on the total i p norm can then be obtained from the p-triangle 
inequality, as follows: 

EEEi^w^E E EEk^^wi p - ( 24 ) 

j m n j QeQj m n 

When p > 1, then the regular triangle inequality will be invoked instead, for instance as in 

(eEEi^^iA <E E (eEk^^>i p ) ' ■ ( 25 ) 

y j m n J j QeQj V m n / 

The rationale for introducing a partitioning into dyadic squares is the technical fact that wave 
atoms are not built compactly supported in space. The windows wq allow to cleanly separate 
different regions of the parameter patch in which the kernel K oscillates with different local wave 
vectors. 

Note also that the dyadic partitioning is a mathematical tool for the proof of Theorem 1, and 
is not part of the construction of the wave atom transform. 

2.3 Geometry of stationary phase points 

Lemma 1 identifies the argument of the Hankel function as a phase. For the kernel K, this phase 
is k(p(s,t) where 4>(s,t) = ||x(s) — x(i)||, and generates typical oscillations as long as V</> has large 
magnitude. In this section, we argue that the locus of near-critical (or near-stationary) points of 
4> necessarily has small measure. The following lemma makes this heuristic precise in terms of the 
scale defect j' . 

Lemma 4. Let (f)(s,t) = ||x(s) — x(i)|| for s,t in some I a . For j' > 0, let 

Kj(f) = {(91,92) : ||V0(2-^ 1 ,2-^ 2 )|| OO < 2"'"}. 
Then there exists C > such that the cardinality of /Cj(j') obeys 

\1Cj(j')\ < C2 j+(j ~ j ">+ 
where (x)+ = x if x > 0, and zero otherwise. 

Proof. Let r = (x(s) — x(t))/||x(s) — x(t) ||; the gradient of the phase is V^>(s, t) = (x(s) -r, — x(t) -r). 
The condition || V<f)(s, f)||oo < 2~ j ' , i.e. 

|x(s) • r| < 2~ j ' and |x(i) • r| < 2~ j> , 

is for large j' an almost-perpendicularity condition between tangent vectors to the curve VL a and 
the chord joigning x(s) and x(t). 

Now fix s = 2~3qi, and let n(s) be either normal vector to tt a at x(s). Let 9 be the angle between 
r and n(s), such that |x(s) • r| = ||x(s)|| |sin#|. Since the parametrization is nondegenerate, the 
first condition |x(s) • r| < 2~i implies 9 < C2~i for some adequately large C > 0. 
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Consider therefore a cone T with apex at x(s), axis n(s), and opening 9 < C2~ 3 ' . The second 
condition |x(i) ■ r| < 2~ 3 is satisfied only if the curve £l a intersects a chord inside the cone at a 
near-right angle, and as a consequence, every chord inside the cone at a near-right angle, differing 
from 7r/2 by a 0{2r 3 ). Because Q a has finite length, bounded curvature, and obeys the geometric 
regularity property (9), there can only be a finite number of such intersections. The total length 
of £l a n r is therefore a 0{2~ 3 ). 

Since the points x(2~ J </2) are a distance C2~ 3 apart from each other, there are at most 
0(max(l, 2 3 ~ 3 )) points indexed by q 2 that obey the two almost-orthogonality conditions, which 
can be written as 0(2^~ 3 )+). Since q\ takes on 0(2 3 ) values, the total number of couples (#1,(72) 
obeying the conditions is 0{2? ' 2^~ 3 " >+ ). 

□ 

We will also need the observation that near-stationary-phase points can only occur far away 
from the diagonal. 

Lemma 5. As before, let d(s,t) = \e 2ms — e 2nlt \. There exist two constants C\,C 2 > such that, 
if d(s, t) < C\, then 

||V#M)||oo > C 2 . 

Proof. As previously, 

||V0(s,t)||oo =min(|x(s) -r|,|x(t) • r|), 

and we write |x(s) • r| as ||x(s)|| | cos(0 s )|, where 6 S is the angle between the chord (x(s),x(i)) and 
the tangent vector x(s). This angle obeys \9 S \ < d(s,t)\, hence the cosine factor is greater than 1/2 
as long as d(s,t) < C\ for some adequate C\. The factor ||x(s)|| is also bounded away from zero 
by regularity of the parametrization. The same argument can be made for |x(t) ■ r|. □ 

2.4 Nondiagonal kernel fragments: decay of individual coefficients 

Within nondiagonal squares, d(s,t) > 1/k and k(p{s,t) > 1, so Lemma 1 asserts that K can be 
written as 

K(s,t) = e ik<t> ^a(k(t)(s,t),s,t), 
where <fi(s,t) = ||x(s) — x(i)|| and the dependence of a on k is mild in comparison to that of e* fc< ^; 

l|l oW (M),M)l<c„ 7 J=«*,tr». 

The presence of additional arguments s and t is needed to account for factors such as the Jacobian 
||x'(i)||; all the derivatives of these factors are 0(1) by assumption. Therefore, the chain rule yields 

a(k<t>(8,t),s,t)\<C a _^ <KM)~ H , #M)<1- (26) 



Fix j > and Qg Qj. We seek a good bound on 



<tfQ,¥W> = / w Q ( S ,t)a(k^s,t),s,t)e ik ^e' i;: —-" 
J3Q 



x 2 3 ipij m \{2 3 s — n\,2 3 t — n 2 ) ds dt, 

where <P{j, m ) has been introduced in equation (3). Without loss of generality, we perform a trans- 
lation to choose the coordinates s and t such that xq = and wq(s,£) = w(2 3 s,2 3 t). 
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A first bound estimating the decay in n can be obtained by using 1) the almost-exponential 
decay (4) for <£>(j im )) 2) the estimate H^qIIl 1 ^ 2~ 2 i that follows from |3Q| < 2 _2jl , and 3) an L°° 
bound for the rest of the integrand, disregarding the oscillations. The result is 

\(K Q ,ip j:m ,n}\<CM2- j sup [(fc<Xs,t)r 1/2 l (l + ||n||)- M , VM>0. 

(s,t)e3Q L J 

The size of the first-order Taylor remainder of k(p(s, t) over 3Q is 0{2 3 ) times smaller that the value 
of k(f)(xQ) itself, so we may evaluate <f> at xq at the expense of a multiplicative constant in the 
estimate. We get 

|(^Q,^,m,n)| < C M 2-- ? '(^( a ;Q))- 1 / 2 (l + ||n||)- M , VM>0. (27) 

Capturing the decay in m, however, requires integrations by parts. Heuristically, the objective 
is to show that the wave atom coefficients decay almost exponentially in m, with a length scale of 1 
in all directions (in units of m), independently of j — at least in the representative case j ~ | log 2 k. 
To this end let us introduce the self-adjoint differential operator 

= I-PA (S:t) -ipk(A4>(s,t)) 
l + P\\kV</>(s,t) - 2im|| 2 ' 

with ^ 

^ = max(2- 2 ^ 2 ,2 2 i)' 

We see that L leaves the exponential exp i[kV (f>(s , t) — 2- ? m • (s,t)] unchanged, hence we introduce 
M copies of L, and integrate by parts in s and t to pass the differentiations to the non-oscillatory 
factors. The scaling parameter (3 has been chosen such that the repeated action of L on the rest of 
the integrand introduces powers of 1/(1 + f3\\kV(f>(s,t) - 2^m|| 2 ), but otherwise only worsens the 
bound by a constant independent of fi = (j, m, n). Indeed, (3 < 2~ 2 \ and 

• the action of each derivative on wq or ^(j.m) produces a factor 2 J balanced by y/J3; 

• the action of each derivative on a produces a factor l/(f)(s, t), which by equation (9) is compa- 
rable to l/d(s, t). Since we are in the presence of nondiagonal squares, l/d(s, t) < min(2 : ', k) < 
2- J . Again, each derivative produces a factor 2- J , which is balanced by \/]3. Note that the lead- 
ing factor 1/ ' \fk^> in the bound (26) is harmless since it is carried through the differentiations. 



It is then tedious but straightforward to combine these observations and conclude that, for all 
M > 0, 

\L M [w(2 j s,2 j t) a(k(/>(s,t), s,t) ^(j, m )(2 j s - ni ,2 j t - n 2 )] | < 

M t/WsJ) (l + /3||feV^,t)-2^m|| 2 )^- 

Since L is a differential operator, the support of the integrand remains 3Q regardless of M, 
hence we still get a factor \3Q\ ~ 2 _2j from the integral over s and t. With the I? normalization 
factor 2 3 coming from equation (3), the overall dependence on scale is 2~ 3 . The resulting bound is 



\{Kq,<Pjw)\<C m 2-* sup (^( S ,t))- 1 / 2 (l + /3||W0( S ,t)-2^m|| 2 )- M 

(s,t)e3Q L 



VM > 0. 

(28) 
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The second factor inside the square brackets can be written as 



fc2-JV^( a ,t)-m 2 \ M 
+ II 2-J/3-V2 11 y ' 

showing that in m-space, it is a fast-decaying bump centered at fc2~ J V</>(s, t) and of characteristic 
width 2-J/3- 1 / 2 . Over the set 3Q, we have the estimate \k2^V(f>(s,t) - k2~^V '4>(x Q )\ = 0{k2~ 2 i). 
The quantity k2~ 2 i is in all cases less than the length scale 2~i f3~ 1 ^ 2 (which is why we could not 
simply have taken (3 = 2~ 2 - ? ), so we may replace V</>(s, t) by V4>(xq) in the expression of the bump, 
at the expense of a multiplicative constant depending only on M. 

We have also seen earlier that (k(f)(s, i)) -1 / 2 can safely be replaced by (/^(xq)) -1 / 2 in the region 
d(s, t) > 1, at the expense of another multiplicative constant. With these observations, we can take 
the geometric mean of (27) and (28) and obtain the central bound 

\(K Ql Vj, m,n)| < C M 2-1 mxQ))- 1 ' 2 (1 + P\\kV<Kx Q ) - 2^m|| 2 r M (1 + ||n||r M , (29) 

for all M > 0. 



2.5 Nondiagonal kernel fragments: £ p summation 

The expression just obtained can be used to show t v summability and verify proper growth as a 
function of k. Only the case p < 1 is interesting and treated in this section. We tackle the different 
sums in the right-hand-side of (24) in the order as written, from right to left. 

• The sum over n is readily seen to contribute a multiplicative constant independent of the 
other parameters j, Q, and m. 

• Consider the sum over m, and pull out the factor 2~i p (k(p(xQ))- p ^ 2 . (We will not worry 
about this factor until we treat the sum over Q.) The range of values for m is an annulus 
1 1 m I |oo x 2 J , so we can compare the sum over m to the integral 



/ (l + /3||A;V</»(x Q )-2^x|| 2 )- M fdx, 



where Cj = {x G R 2 : d2 j < < C 2 2 j } for some C\,C 2 > 0. In what follows, take M 

sufficiently large so that, say, Mp > 5. Two cases need to be considered, corresponding to 
2 2j < k (large scales), and 2 2j > k (small scales). 

— If 2 2j < k, then j3 = 2 2j, fc -2 . It will be sufficient to consider only the upper bound for 
1 1 x I |oo, whence we have the bound 



/ 



-Mp 

I 1 + || f -KW V ^Q) ~ A. 2 \ 

J <C2i ' ■->-->" ■ 

With Lemma 4 in mind, we introduce the scale defect j'g as the unique integer such that 

V^<||V<K* Q )||<2- J Q. 

The integrand is a bump that essentially lies outside of the region of integration as soon 
as k2- {j+j 'Q ) > 7?. 



15 



More precisely, observe that 

-Mp 



\\oo<CV \ 2 3k J \ ) 



xiUxll^C^' 

hence the integral is bounded by a first expression, 

-2Mp 



c 



2 2j (l + 2 j (2- J Q - Ck~ l 2 2j ) + ^ P . (30) 



A second bound can be obtained by letting x' = x — 2 J kV(p(xQ) and extending the 
region of integration to the complement of a square in x', of the form 



\Aoo>s jQ =(±k2-v + W-cvY 



If SjQ < k2~ 2 i , we might as well put it to zero and obtain the bound C{k 2 i) 2 for the 
integral. If SjQ > k2~ 2 i, the integrand can be made homogeneous in x and the integral 
bounded by 



Now uniformly over SjQ, the resulting bound is 

C{k2- 21 ) 2 (l + 2 3 {2- 3 Q - fe- 1 2 2 - ? )+J . (31) 

The minimum of (30) and (31) is 

C min (2 2j , k 2 2~ 4j ) (l + 2 j (2- j Q - fc" 1 2 :y ) + )~ 2MP+2 . (32) 



— If k < 2 2j , then [5 = 2 2j . This time we will only consider the lower bound for ||x| 
and write 

(1 + \\2- j kV(/)(x Q ) - x|| 2 )- Mp dx. 



L 



x\\oo>C2i 

Since ||V</>(xq)|| is a 0(1) and 2~ik < 2 J , there exists a value j* < \ log 2 k + C such 
that for all j > j*, the center of the bump is inside the square ||x||oo < C2- ? (the 
constant C changes from expression to expression.) When this occurs, we can let x' = 
x — 2~i kV(fi(xQ) as before and consider the integral outside of a smaller square and 
bound 

(l + ||x|| 2 )- Mp dx<2- 2 ^ Mp - 2 \ j>j*. (33) 



L 



For the few values of j such that \ log 2 k < j < j*, we recover the previous estimate, 
namely C (/c2 _2 - J ) 2 , which is a 0(1). 

Consider now the sum over Q, and recall that the bounds just obtained need to be multiplied 
by 2~i p (k(j)(xQ))~' p l 2 . Pull out the factor 2~~ JP one more time. Again, we need to separately 
consider 2 2j < k (large scales) and 2 2j > k (small scales). For small scales, the bound (33) is 
uniform in Q, hence the sum over Q £ Qj simply contributes a factor 2 2j . 
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For large scales, the strategy is to split the sum over Q into a near- field contribution, for 
which d(s, t) < C\ in the sense of Lemma 5, and a far-field contribution, for which stationary 
phase points must be handled adequately. The terms in the far-field sum is then further 
broken down into groups corresponding to a same value of the scale defect, which Lemma 4 
helps identify. Schematically, 



E = E + E 

<9eQi Qe near-field f>o 



E 

Q-. scale defect=j' 



Consider the two regions separately. 



Near-field. In this region, 4>(xq) may be as small as 1/k, hence we estimate (k(f>(xQ)) p l 2 < 
1. By Lemma 5, the scale defect j'g is bounded by a constant, hence the bound (32) 
becomes 

C min (2 2 '',fc 2 2- 4 '') (l + 2^(1 - AT^')^" 2 ^ 2 . (34) 

We claim that this quantity is always less than a constant independent of j and k. 
Indeed, if j is so large that 2 2j k' 1 > 1/2, then k 2 2~ 4j < 4, and it suffices to use the 
trivial minoration 1 + 2- ? (l — k~ 1 2 2] ) + > 1. If on the other hand 2 2j fc~ 1 < 1/2, then we 
have 1 + 2 j (l - k~ 1 2 2j ) + > \2 3 , which implies that (34) is bounded by 

C min (2 2 ^, A: 2 2- 4 ^) 2~ 2 '^ M ^ < C 2' 2 ^' l+M ^ < 1, 

because we chose Mp > 5. The sum over Q then contributes a factor proportional to 
the number of nondiagonal squares, i.e., 

Far-field. The leading factor (k4>{xQ))~ p l 2 now contributes a factor k~ p l 2 , since 4>{ x q) ^ 
C in the far field. For the sum over Q, we use equation (32) one more time and write 



C min (2 2 ^', fc 2 2-^') ^ ( X + 2 ^ 2 ^ " k ^ 22l )+) 
Qe far-field 



-Mp+2 



For each Q, find the closest integer j' < j to jg. As long as j' < j, Lemma 4 asserts 

that the number of terms comparable to ^1 + 2 J (2 -J ' — k~ 1 2 2 i) + ^ P is aO(2 2 ^'). 
The endpoint j' = j receives the contribution of arbitrary large jg, meaning terms that 
can be as large as a O(l); however by Lemma 4 there can only be 0{2 ] ) such terms. 
After indexing terms by j' in place of Q, we get the bound 



Cmin {2 2j ,k 2 T 



4j\ 



23 + ^ 2 2j ~ j ' (l + 2 j (2^" - A; _1 2 2j ') + ) MP+2 

-C<j'<j 



It is easy to see that the summand peaks for f near jo = min(j, — 2j'+log 2 k); it decreases 
geometrically for j < jo because of the factor in brackets, and decreases geometrically 
for j > jo because of the factor 2 2 - ? ~- 7 . The result is a bound 

C min (2 2j , k 2 2- Aj ) max^'fc -1 , 2 j ) = C mm(k, 2 3j ). 
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• What remains after gathering the various bounds is a constant Cm,p times 



2~ jp 2 2j (near-field) 

i<|log 2 k+C 

+ k- p/2 2- jp mm(k,2 3j ) (far-field) 

i<|log 2 k+C 

+ J2 2 2j 2~ jp 2~ 2 ^ Mp - 2 \ (small scales) 

j>|log 2 k+C 



The near-field contribution sums up to C p k 1 ' p l 2 as soon as p < 2. The far field contribution 
is bounded by 



2 i(3-P) + k Y 2- jp 

j>l l°g2 k 



< Cpk'P^k 1 -^ 3 < c p k 1 - p/2 . 

With the choice Mp > 5, the contribution of "small scales" is negligible in contrast to the 
first two terms. 0{k 1 ~ p / 2 ) is the desired growth rate in k, compatible with equation (19). 
This concludes the part of the proof related to nondiagonal squares. 



2.6 Diagonal kernel fragments: decay of individual coefficients 

It is now assumed that the dyadic square Q overlaps with the diagonal strip S = {(s, t) : d(s, t) < 
1/k}. Because of the singularity of the kernel at s = t, the integrations by parts cannot proceed as 
before. Inside S, the smoothness of the kernel is governed by the case kx < 1 in Lemma 1. 
Further complications arise, depending on the value of the scale j. 

• At scales j < log 2 k + C, a square Q intersecting with S is not entirely contained in S; in 
fact, a large portion of it lies in the nondiagonal portion d(s,t) > 1/k. We call this portion 
(two triangles) the off-strip contribution. There are 0(2- ? ) such triangles. 

• At scales j > log 2 k+C, some squares may be contained inside the strip S without intersecting 
the diagonal s = t. These squares make up the regular on-strip contribution; there are 
0{2 2 ^k^ 1 ) such squares. 

• The remaining 0(2 J ) squares or portions thereof, overlapping with the diagonal s = t, make 
up the singular on-strip contribution. 

In order to smoothly cut off the strip S from dyadic squares, introduce a = s — t (defined 
modulo 1 in [0, 1]) and r = s + t. By symmetry of the problem, one can consider the triangle 
{(s,t) e Q : s > t} and still call it Q, without loss of generality. We can therefore focus on a > 0. 
If we properly select the coset relative to the modulo operation, we can also take (<r, r) to smoothly 
parametrize the triangle Q. With these choices, the diagonal strip is S = {(a, r) : a < 1/k}. 
Consider now a smooth indicator p(ka) where p is a C°° positive function obeying 

, , f 1 if < x < 1; 
P{X) = \ ifx>2. 
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Multiplying the integrand in (41) by p(ko~) gives the on-strip contribution (regular and singular); 
multiplying it by 1 — p(ko~) gives the off-strip contribution. These cases are treated separately. 

• Off- strip contribution. Over the off-strip region we have ka > 1 hence k(f>(s, t) > 1 by equation 
(9), so that the case kx > 1 of Lemma 1 applies there. The analysis of the coefficient decay 
in n is the same as in the previous section, so we omit it here. 

As far as analysis of the decay in m, we are back in the setting of the analysis of Section 
2.4, except for the factor 1 — p{ka) that prevents the same scheme of integrations by parts in 
a. Each derivative of p{ka) would produce an unacceptably large factor k. The smoothness 
in r is however unaffected, which permits to carry over the analysis of Section 2.4 with 
integrations by parts in r only. The direction of increasing r in the m plane is e T = (1,1), 
to which corresponds the decomposition 2 m ■ (s,t) = (mi + ttt^t + (mi — 777,2)0". Since 
V</>(0) = lim cr ^o+ ^4> points in the direction of a, we have e r • V0(O) = 0, and we can use 



L = 



1 + /32 2 J («n±im) 2 



to generate the repeated integrations by parts, where (3 is the same as previously. Notice 
that the amplitude is uniformly bounded, since kx > 1 in the off-strip region. The result is 
a bound that involves mi + 7772 only. With the contribution of the decay in n, the off-strip 
coefficient estimate is 

I (Kq, (1 - p(ka))(p j:m , n )\ < C M 2~ j (1 + /32 2j (mi + ™ 2 ) 2 y M {\ + ||n||)- M , (35) 
for all M > 0, and only for scales obeying 2 2j < k. 

On-strip contribution: amplitude estimate. We now take ka < 1. The case kx < 1 of Lemma 
1 allows to write 

K(s,t) =e ik<t>{s ' t) a(k ( f ) (s,t),s,t), 
where now the amplitude's smoothness is 

d n a 

\^(kcP, S ,t)\<C n <p-\ <f><l. (36) 

The partial derivatives of a with respect to the arguments s and t are O(l) and well within 
the above bound as long as </>(s,t) < 1. To compute the total derivatives with respect to s 
and t, however, it is necessary to contrast smoothness along and across the oscillations, by 
means of the coordinates a and r. The value of (f>(s, t) is comparable to the circle distance 
d(a, 0) , namely 

Dd(a,0)<cf>(^,^^j <Dd(a,0). 

The first inequality is exactly equation (9), the last inequality follows from a Taylor expansion. 
Since we only consider a > 0, we write this property as <j) x °~- A careful analysis of Taylor 
remainders shows that the same estimate estimate is true for the r derivatives, 

d n 4> f a + t t - a\ 

<C n a, o-^O, (37) 



< C n , + 0. (38) 



'dr 1 

while the a derivatives do not yield any gain: 

. d n (f> ( a + r t — a 
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The action of the successive r derivatives on a through its 4> dependence can be understood 
from the higher-order analogue of the chain rule, known as the combinatorial Faa di Bruno 
formula: 




Tren 



dr\ B \ " 



In this formula, II is the set of all partitions tt of {1, . . . , n}; | vr | denotes the number of blocks 
in the partition tt; these blocks are indexed as B G tt; and \B\ denotes the size of the block B. 
Since there are \tt\ factors in the product over B, equation (37) reveals that the derivatives 
of (f> yield a factor C n a^. On the other hand, by equation (36), each ^differentiation of a 
introduces an inverse power of a. The order of the derivative is \tt\, for a contribution of cr _ l 7r l 
that exactly cancels the coming from the derivatives of (f>. 

This analysis only concerns the dependence of a on r via (j). It is easy to apply the multi- 
variable chain rule to see that the depedence of a on r via its second and third arguments (s 
and t) does not change the conclusion that any number n of r derivatives, n > 1, keep the 
amplitude bounded, with bound independent of j and k (but not n, of course): 

'^K^^)'^^) I£C ~ »^W0. (39, 

There are no factors to gain in the a derivatives of the phase, hence the same analysis yields 
*>af..f* + T r-a\ „ + TT -<A.^_ B> ^ < x , n > W . (40) 



'da 1 

We are now equipped to study the coefficient 

(Kq, <p jt m , n ) = / w{Vs, 2^)e^ s '*)a(^( S , t), s, t) e ~ l2Jm < s ^ 
•J3Q 

x 2 j ^(j,m)(2 j s - m, 2 j t - n 2 ) (is eft. (41) 

• Regular on-strip contribution 

For regular on-strip squares, i.e., those squares at very small scales j > log 2 k + C that 
intersect with the strip S but not with the diagonal s = t, the decay in m and n is obtained 
by a simple argument of integration by parts. Introduce copies of 

and integrate by parts in (41). Each derivative in a acting on the amplitude a(k(j),s,t) 
produces a factor x cr _1 . Since Q does not intersect with the diagonal, a > hence 
a l ^5 This factor is balanced by the choice of scaling in the expression of L 2 . A fortiori, 
the derivatives in r are governed by a stronger estimate and are therefore under control. 
The derivative in r or a acting on <^(j )iri)n ) do not compromise its super-algebraic decay, 
hence we gather the same decay in n as previously. One complication is however the possible 
logarithmic growth near a = of the amplitude a when it is not differentiated. Consider the 
intermediate bound 

\(K Q ,p(ka)ip jiia , n )\ < C M V (l + ||m|| 2 )^ M (l + ||n||r M [ |a(fc^(s, t), s, t)\ dsdt, 

J3Q 
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for all M > 0. If the amplitude were bounded, then the integral would produce a factor 2~ 2 i 
like in the nondiagonal case. Instead, we claim that the integral factor is bounded by 2~ J A; -1 . 
In order to see this, consider the bound 

\a(k<f>,s,t)\ <C(l+|log(fc0(a,t))|), 

from Lemma 1. Since log is increasing and <\> x a, there exist C\, C 2 > such that 

l0g(Clfc<7) < log(fc0(s,t)) < log(C 2 fc<7), 

hence 

|log(fc0(a,*))| < C+ I log(fccr)|, for some C > 0. 

This bound does not depend on r, and since (s, i) G 3Q, r ranges over a set of length 0(2 -J '). 
We therefore obtain the bound 

1 

f \a(k<t>(s,t),s,t)\dsdt<C2-i x (C + | log(fe<r)|) da < C^AT 1 . 

The final estimate for the regular on-strip contribution is 

\(K Q ,p(ka)<p jim:n )\ < C M k- 1 (1 + ||m|| 2 )~ M (1 + ||n||)- M , (regular on-strip) (42) 

• Singular on-strip contribution 

Let us now consider a dyadic square Q that intersects with the diagonal s = t. Heuristically, 
we cannot expect that the decay length scale of the wave atom coefficients be independent 
of j in all directions in m: because Q overlaps with, or is close to the diagonal, the decay in 
the direction mi — mi is much slower than the decay in the direction m\ + m^. However, the 
number of diagonal squares is small enough to restore the overall balance at the level of the 
£ p summability criterion. 

To quantify the decay in the m\ + m^ direction, introduce the self-adjoint operator 

I-2- 2 ^ 

t Or A 

l+(^) 2 

It leaves the exponential e~ l2J ( m i s + m 2t) invariant. After integrating by parts, the action of 
/ — 2~ 2 - ? ^ 3 leaves the bound on the rest of the integrand unchanged, because 

(i) w(2is, 2H) and v ? 0,m)(2"' s — ni,2H — ni) produce a factor 2 J when differentiated; 

(ii) differentiating (/ , (j,m)(2-' s — n\,2H — 77-2) does not compromise its decay in n; and 

(iii) a(k(f), s, t) has a logarithmic singularity, and otherwise becomes uniformly bounded when 
differentiated in r, as we have seen. (The presence of the scaling 2~ 2 i in L\ is not even 
needed here.) 

The integral J 3 q \a(k<p(s, t), s, t)\ dsdt for the amplitude can be bounded by C 2~- ? max(2 _J ', k^ 1 ) 
as we arued for the regular on-strip squares (here 3Q is not necessarily contained in S.) The 
result is a bound 

\{K Q ,p(ko-)<p jtm , n )\ < C M max(2^,^ 1 ) (l + (mi + m 2 ) 2 )~ M (l + ||n||)- M , (singular on-s 

(43) 
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for all M > 0. 

Finally, the decay in mi — m,2 for those (singular, on-strip) squares that intersect the diagonal 
cannot proceed as previously An analysis of coefficients taken individually would be far 
from sharp, e.g., would not even reproduce £2 summability. The proper reasoning involves a 
collective bound on the £2 norm of all the wave atom coefficients at a given scale j > 0, which 
correspond to squares Q that intersect with the diagonal. This reasoning is explained in the 
next section, and gives the bound 

^ \{Kq, p{ko)(pj )m ^)\ 2 < C j 2 2~ 3j , Q G Qj and Q intersects the diagonal. (44) 

m,n 

The study of £ p summability from all these estimates is then treated in Section 2.8. 
2.7 Diagonal kernel fragments: collective decay properties 

The strategy for obtaining (44) is to compare wave atom coefficients to wavelet coefficients, scale 
by scale. Estimating individual wavelet coefficients is a much tighter way to capture the sparsity 
of a log singularity than directly through wave atoms. (Wavelets however are not well- adapted for 
the overwhelming majority of dyadic squares that correspond to C°° oscillations.) 

Consider two-dimensional compactly supported Daubechies wavelets with one dilation index j, 
built on the principle of multiresolution analysis [21]. They are denoted as 

r f , n (s,t) = V'i; £ (2i's - ni ,2?t - n 2 ), 

where e = 1, 2, 3 indexes the type of the wavelet (HH, HL or LL). The easiest way to define Meyer 
wavelets in a square it to periodize them at the edges. 

Fix j < 0, consider a function f(s, t) defined in [0, l] 2 , and consider its wave atom coefficients at 
a scale j. By Plancherel for wave atoms, there exists an annulus Aj = {(^1,^2) : C\2 2 i < ||£||oo < 
C 2 2 2j } such that 

Ei(/.to)P< / \f(0\ 2 d£. 

m,n J 

By Plancherel for wavelets and the properties of Daubechies wavelets [21], there exists jo such that 
this L 2 energy is for the most part accounted for by the wavelet coefficients at scales 2 j — jo < j' < 
2j + jo , i.e., 

\f(0\ 2 dt;< £ El(/-#,n)l 2 . (45) 

j'e[2j-j ,2j+j ] e,n 

The last two equations show that, collectively in an £2 sense, wave atom coefficients at scale j can 
be controlled by wavelet coefficients at scales neighboring 2j. 

The relevant range of scales for this analysis is j > ^ log 2 k. The on-strip region has length 
0(2 _J ) and width 0(min(/c _1 , 2~-?)). Since each wavelet is supported in a square of size ~ 2~ 2j, -by- 
2~ 2j , the number of wavelets that intersect the strip is 0{2? x2 2 i min(A; -1 , 2 -J ')) = 0(min(2 3j A; -1 , 2 J ')) 
Among those, only 0{2 ] ) correspond to wavelets intersecting with the diagonal s = t. The bound 
on wavelet coefficient depends on their location with respect to the diagonal: 

• Non-diagonal wavelets. The wavelet's wave number is ~ 2-?" ~ 2 2 ^' and soon becomes much 
larger than the local wave number ~ k of the oscillations of the kernel, hence a fast decay in 
j' — > 00. More precisely, fix Q 6 Qj; the coefficient of interest is 

{K Q p(ko-),ip £ jln ) = [ a(k^(s,t),s,t)e ik ^ t) w(2 j s,2 j t)2 j ''tp £ (2 j 's-n 1 ,2 j 't-n2)dsdt. 
./supp n 
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Because the wavelet has at least one vanishing moment, one may write it either as 
iP £ (2 j 's - m,2 j 't - n 2 ) = 2- j,< ^-(2 j 's - ni ,2 j 't - n 2 ), e = HL or HH, 

or as 

tP £ (2 j 's - n u 2 j 't - n 2 ) = 2^'^(2 J "s - n u 2 j 't - n 2 ), e = LH, 
where ip £ has the same support as ip £ . After integrating by parts in s or t, we can 

- use the bounds (39) and (40) on the amplitude; 

- use the bound V (S)t) e^ = O(k); 

- use V {Sit) w(2is,2H) = 0(2i); 

- use |supp ^ £ , n | < 2- 2 ^' 

to conclude that the coefficient obeys 

\{K Q p(ka),^ l>n )\ < C M 2~i' 2~i' (max(k, 2>') + ^) , 

where cf)" 1 is a notation for the supremum of <p _1 over the support of the wavelet. If we 
index by the integer q > 1 the distance between the center of the support of the wavelet 
to the diagonal, as \[2q2~i , then (j)' 1 x q~ l 2^ . The bound above becomes Cm 2~ J {q" 1 + 
2-J" max(ife,2^)). 

The sum in the right-hand- side of (45) is then estimated as follows. As we saw earlier the 
length of the strip S n Q is 0(2 _JI ), and its width is min(l/A;, 2~ J ). Since the translation 
step of wavelets is 2~3 ~2"^, the translation index n takes on 2 J x 2 2 i min(l/fc, 2 values. 
Hence 

£ ^K^Qp(fea),#,n)| 2 <CEl 2 ^^ 1+2 ^ fc )| 2 
j'e[2j-j ,2j+jo] £ > n n 



< C2 J ' ^(2~ 4 V 2 ) + 2 3j ' min Q, 2"^ x 2- 4j '(2~ 2j ' max(fe, 2 J ')) 2 



9>1 

Q ■ l 

< C2^ 3 . (because j > - log 2 k + C) 

• Diagonal wavelets. For wavelets intersecting the diagonal, it will not be necessary to quantify 
cancellations. By Lemma 1, 

\(K Q p(ka),^ £ , n )\<C [ {l + \\og{k<j)\)2i' \^ £ {V' s-rn^'t-nz^dsdt. (46) 

./supp n 

Without loss of generality we can consider ka < 1/2 and write 

| log 2 (fcff)| = - log 2 (/c(j) = - log 2 (2 J " a) - log 2 k + f < - log 2 (2 J "a) + /. 

Since log is integrable near the origin, and |supp ipj/ n \ ^5 2^ 2 ^ , the contribution due to 
— log 2 (2 J, 'o") is a 0(2~i'). The contribution of the lone f, on the other hand, is a 0(j'2~ 3 '). 

There are 0(2 J ) diagonal wavelets, each with coefficients 0{j'2~ 3 ') = 0(j2~ 2 i), hence the 
sum of their squares in the range 2j — jo < j' < 2j + jo is 0(2? x (j2 _2j ) 2 ) = 0(j 2 2 -3j '). 

As j — > oo, the contribution of diagonal wavelets manifestly dominates that of non diagonal 
wavelets, and we have shown that the resulting estimate is (44). 
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2.8 Diagonal kernel fragments: £ p summation 

Let us conclude by calculating the growth of Ylj Y^QeQj S m n K-^Qj c / 3 i,m,n)| p in the parameter k, 
for those dyadic squares that intersect the strip ka < 1 . We start by letting p < 1 . 
Consider each contribution separately. 

Off-strip contribution. Recall that j < \ log 2 k + C in this case. The reasoning entirely 
parallels that of the previous section and we encourage the reader to focus on the discrepancies. 
First, the sum over n yields a harmless constant factor. Second, use (35) and drag the factor 
2~ JP out of the sums over m and Q; the former sum is then comparable to the integral 



1= [ {l + \/3 1 / 2 2i(x 1 + x 2 )\ 2 )- M rdx 1 dx 2 , 

JC-i 



with Cj an annulus of inner and outer radii proportional to 2- J . Over this domain, the 
integrand concentrates near the union of two "ridges" of length ~ 2 J and width ~ (3~ 1 l 2 2~i , 
oriented along the anti-diagonal x\ = — x-i. Note that (3 = 2 2j /c~ 2 . The integral / is therefore 
bounded by a constant times 2 J x (/3 _1 / 2 2 _J ) = 2~- ? /c. Third, the sum over Q G Qj that 
intersect with the strip yields a factor 2 3 , proportional to the number of diagonal dyadic 
squares at scale j. The remaining sum is bounded by 

C p 2 j (2~ j k)2- jp < Cpk 1 '^ 2 . 

i<|io g2 k+c 

This is the desired growth rate in k. 

Regular on-strip contribution. Here, j > ^log2 k -\- C , therefore HmHoo > C2i >CVk. The 
factor (1 + 1| rci||oo) ^ in equation (42) therefore yields a negative power k ^/ 2 for all M > 0, 
i.e., /c~°°, that sums up to a negligible contribution. 

Singular on-strip contribution. As previously, two scale regimes should be considered. When 
j < 2log2&, we can use the bound (43). The sum over n is harmless; the sum over m 
produces a factor 2 J since there are significant 0(2- 7 ) values of m on the ridge \mi +7712I < C 
at scale j; the sum over Q produces another factor 2 J since there are 0(2- ? ) diagonal dyadic 
squares a scale j. The resulting sum over j is then bounded by 

C p Yl 2 2j 2~ jp < Cpk 1 ^ 2 , 
i<\ i°g2 k 

which is again the desired growth rate. 

If now j > I log 2 k, we need to invoke the collective decay estimate (44. Fix j and Q £ Qj 
a singular dyadic square. By equation (43), values of mi + 771,2 significantly different from 
zero will give rise to negligible coefficients that sum up to a o(k). More precisely, let 5 > 
be arbitrarily small. Then the wave atom coefficients in the region \mi + 771-2! — C2 s i 
decay sufficiently fast (take M S> 1/(5) that their total contribution is a 0{k~°°) in £ p . The 
significant coefficients at scale j are, again, on a ridge of length 0(2'') and width 0(2 s i), for 
a combined total of iV = 0(2^ 1+<5 )) significant coefficients. 

We can now relate the £2 norm estimate (44) to an l v estimate, < p < 2, by means of the 
Holder inequality 

/ \1/P / \l/2 

[J2\{KQ P (ka), Vj , m , n )\P < \ ^2\(K Q p(ka),^ mtn )\ 2 x N*~*, (47) 
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where N = 0{2^ 1+ ^). After simplification, the right-hand side is bounded by Cs, p 2 + 1 /p+ s ') 
where 5' is another arbitrarily small number, namely 5' = — |). This quantity still needs 
to be summed over Q — there are 0(2 J ) such squares — and then over j; but the summation 
method will depend how p compares to 1, and accordingly, which of equations (24) or (25) 
should be used. 

If p < 1, then (24) should be used, and we obtain 

ll^lli(,)<^ E 2^2^ 2 +^')) P = C 5 , p £ 

which always diverges. However if p > 1, then (25) implies 

i>5 log 2 k 

For any 1 < p < 2, the above series is convergent if for instance we choose & = \{^- — ^j, 
and the result is a O(l), independent of k. This part of the singular on-strip contribution 
falls into the second category identified at the beginning of the proof, namely equation (20). 
This concludes the proof in the case when K is the single- layer potential Go- 



2.9 Analysis of the double-layer potential 

The proof of the sparsity result for the double-layer kernel G\ defined in equation (11) is a simple 
modification of that for the single-layer kernel Go- 

• Nondiagonal part. The smoothness bound for Hankel functions in Lemma 1 exhibits the same 
rate for all n > in the case when x > 1/k. The other factors, functions of s and t which 
accompany the Hankel factor in formula (11), have no bearing on the sparsity analysis since 
they are smooth and do not depend on k but for the leading factor ik/4. As a consequence, 
G\ can still be written as a product 

Gi(s,t) = a(fc0(s,t),s,t)e ifc<?!,(s ' i) , (nondiagonal part, \s-t\>l/k), 

where the amplitude obeys the same estimate as previously, but for a factor k: 

\-^a(kcj)(s,t),s,t)\ <C n k—^====4>(s,t)~ n . 

With a number of nonstandard wave atom coefficients |A| = O{ke" l / co ), one could form an 
approximation of the nondiagonal part of Go with error e. In the case of G±, the same number 
of terms results in an error that we can only bound by ke. Thus, to make the error less than 
a specified e, the number of terms needs to be 0(fc 1+1 / 00 e~ 1 / 00 ). This justifies the form of 
the first term in (14). 

• Diagonal part. For x < 1/k the smoothness estimate in Lemma 1 is worse for than for 

but as is well-known, the dot product 

x(s) — xfi 
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is small near the diagonal and more than compensates for the growth of there. The 
precise version of this heuristic is a decomposition 



G\(s, t) = a(k(p(s, t), s, t) e 



ik<j>(s,t) 



(diagonal strip, |s — t\ < 1/k), 



where we claim that the amplitude obeys the same estimates as those for Go in variables 
a = s — t, t = s + t, namely 



a + r t — a\ a + T T — a 



d n a 
At™ 



(Pa 

a + T t — a\ a + T T — a 



< G n , 



< c n < 



<1, 



(48) 



(49) 



for a < 1/k, and this time for every n > including zero. (Total derivatives in a are defined 
by keeping r fixed, and vice- versa.) 

Let us prove (48). As previously, put r = ||x(sjlx(t)|| • Observe that n x ( t )||x(t)|| = (x(t))- 1 • r. 
Derivatives of (x(t))^ ■ r in a and r are treated by the following lemma, proved in the 
Appendix. 



Lemma 6. For all n > 0, 



*dr" 
d n - 
'do" . 



(x(t)^.r 



I < C n a, 
<C n a l - n . 



(50) 
(51) 



The chain rule and Faa di Bruno formula can then be invoked as previously, with the combined 
knowledge of (50), (51), the growth of from Lemma 1, i.e., 



d" 



'd(f) r ' 



^ 1) (fc(/>(s,t))e" 



-ik<t>{s,i) 



<C m — a 
ka 



as well as equations (37) and (38) on the growth of the derivatives of (f>. It is straightforward 
to see that (48) is satisfied; for instance, the factor a from (x(t))- 1 • r and the leading k in 
the expression of G\ cancel out the l/(ka) in the formula for the derivatives of H^. The 
rest of the argument involving the Faa di Bruno formula is the same as previously. Equation 
(49) follows from the same reasoning, and the observation that no a factor is gained upon 
differentiating <p in a. 

Since (48) and (49) are at least as good as what they were in the case of Go, the rest of 
the argument can proceed as previously with the same results; the "off-strip" and "regular 
on-strip" contributions, for instance, are unchanged from the Go scenario. The "singular 
on-strip" contribution however, corresponding to dyadic squares that intersect the diagonal, 
ought to be revisited since G± has a much milder singularity than Go near the diagonal. 

The estimate of fast decay in \m\ + 7712 1 and ||n||, namely (43), is a fortiori still valid. It 
appears, however, that the collective bound (44) at scale j can be improved to 



^\(K X f ag p(ka),^ 



< C 2- 6j k 2 , 



(52) 
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where xf '^ ( s > t) refers to the wq(s, t) over the squares Q G Qj at scale j for which the 

support of wq intersects the diagonal. The presence of an aggregation of windows Xj^^ ( s > i) 
is important here, as the study of coefficients corresponding to individual windows wq would 

not give a sharp bound. Whether xf^ & ^ i s > t) or p{ka) effectively determines the cutoff 
depends on the relative values of j and log 2 k. 

Again, via a Plancherel argument, the scale-by-scale bound (52) can be proved by passing to 
a system of Daubechies wavelets. We have 

{K X f^p(ka),r r , n ) = I xf ag (s,t)p(ka) (x(i))^ • r (53) 

Jsupp^, in 4 

x2 j ' ij £ (2 j 's - m,2 j 't - n 2 ) dsdt, 
where the scale of the wavelet relates to that of the window w as 2j — jo < j' < 2j + jo. 

As previously we will use the vanishing moments of the wavelet to bring out a few 2 _J ' factors. 
This time we will need up to three vanishing moments, i.e., we write the wavelet as 

/ j \ M 

i) £ (2 j ' s-n l ,2 j 't-n 2 )= h- j ' —\ 4> £ {2 j ' s - n u 2 j 't - n 2 ), e = HL or HH, 

or as 

/ , \ M 

^ £ {2 j 's-n 1 ,2 j 't-n 2 ) = h- j '-j 4> £ {2 j ' s - m,2 j 't - n 2 ), e = LH, 

where M < 3, and i\) £ has the same support as ip £ . Before we let these derivatives act on the 
rest of the integrand, multiply and divide by 4>(s,t) = ||x(i) — x(s)|| to get 

kH^ik^s, t)) (x(t)) x • r = \k<f>(8, t)H^\k4>(s, t))} (x(t)) ± - X(t) ~ X(S) 



|x(t)-x( S )||2- 

We will then need the following lemma, which refines the results of equations (37), (38), and 
Lemma 6. It is proved in the Appendix. 

Lemma 7. Let 4>(s,t) = ||x(s) — x(i)||. For every integer m < there exists C m > such 
that, as long as s ^t, 

' ' <Ks,t)\<C m , (54) 



ds 



ds 



(i(t)) ± x(t)-x( a ) 



< C m , (55) 



|x(t)-x( S )PJ 

The same inequalities hold with d/dt derivatives in place of d/ds derivatives. 

Let us first consider the wavelets whose support intersects the diagonal. There, the support 

of the wavelet is sufficiently small that Xj^ a ^( s > t) = p{ka) = 1. One may integrate by parts 
only once in (53), because 

V (i)t) ^(s,t) = (x(s)-r, -x(t) t) 

is discontinuous, since the unit chord r = (x(s) — x(t))/||x(s) — x(t)|| changes sign across 
the diagonal and x(i) / there. The action of either d/ds or d/dt on the integrand after 
integration by parts gives: 
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— a 0(1) contribution for -^(xH[ (x)), by Lemma 2; 

— by the chain rule, a 0(k) contribution for ^{k<p) and ^(k<p), because of Lemma 7; 

— a 0(1) contribution for derivatives of (x(i))- 1 • r/(f>, by Lemma 7. 

The size of the support is 2~ 2 - ? ' , the wavelet comes with an L 2 normalization 2 J ', one factor 2 - - 7 ' 
comes out of the vanishing moment, and \ f — 2j\ <const.; hence diagonal wavelet coefficients 
obey the bound 

\(K X f ag p(fca), V/, n )l < 2 ~ 4j fc- (diagonal wavelets) 

There are 0(2 J ') = 0(2 2j ) such diagonal wavelets overall, hence the sum of squares of these 
coefficients is bounded by 2~ 6: >k 2 , in accordance with equation (52). 

Let us now treat the wavelets that do not intersect the diagonal s = t, and show that the 
same bound is valid. One will now need to integrate by parts three times in s or t to get 
three 2~ 3 factors out, and gather the action of the derivatives on the rest of the integrand as 
follows. 

— The factors p(ka)Xj^ a ^(s, t) are essentially multiplied by max(2 J , k) for each derivative. 

— By Lemma 2, the combination xH^\x) becomes 1/x when differentiated three times in 
x. This is \/k<p when x = k(f>. 

— Derivatives of (p in s and t remain O(l) by Lemma 7, hence each derivative of k<f> yields 
a factor k. 

— By Lemma 7, al derivatives of (x(t))- 1 • [j~^~~^jp remain bounded uniformly in s and 
t. 

The product rule yields many terms but the overall sum is controlled by the behavior of the 
"extreme" terms identified above, hence a factor max(2 3j , k 3 ) + k 3 /k(p under the integral sign. 
Since the wavelet has support well away from the diagonal, we can proceed as previously and 
bound (/> -1 (s,i) by g~ 1 2 J where q is an integer indexing the distance between the diagonal 
and the center of the wavelet. Again, the support of the wavelet has area 0(2- 2 ^") and f is 
comparable to 2j, hence we get a bound 

l(^xf g /'W.^p)l < 2 ~ 8j [max(2 3 - ? ', k 3 ) + k 2 q~ 1 2 2j ] . (nondiagonal wavelets) 

As seen previously, the number of wavelet coefficients is a 0(min(2 2:, /A:, 2 J )) across the diago- 
nal (indexed by q), times 0(2 2 - 7 ) along the diagonal, for a total of 0{vci\n{2 A ^ /k,2^)). Hence 
we have 

£ ^|(K Q p(^),^,„)| 2 

j'e[2j~j (h 2j+j ] s,n 

2 2 ^2- 16 ^V 2 2 4j 
q 

< [2- 12 ik 5 + 2-^'j + [2- 10 ^ 4 ] 

< 2- 6 ^ 2 since j>- log 2 k + C. 

This is the desired decay rate, compatible with equation (52). 



< min (^-,2 3j ^j x 2~ 16 '?max(2 6j ',A; 6 ) 
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We are now left with the task of verifying that (52) implies the correct decay of the l p norm, 
as in equation (21). Let p < 1. Start by using Holder's inequality (47) with N = 0(2^ 2+<5 ))- 
there are 0(2 J (1 + 5)) wave atoms per square Q, and 0(2 J ) squares along the diagonal. We 
get 

J2\{K X f &g p(ka),^ mtIl )\P < C ^(2-3^) 2 ^2+<5)(i/p-i/2)] p = c u> 2 Q-Wrt, 

m.n 

where 8' = S(^ — Finally, the p-triangle inequality asks to sum this bound over j > | log 2 k. 
The sum is convergent provided p > 1/2 and 5' is taken sufficiently small. The result is 

l|tfj£ p(F) < C p k p k^ 2 - 4p+5 '^ = C p k l - p+s \ V<5" > 

After taking the 1/p-th power, we fall exactly into scenario 3 for the £ p summation, i.e., 
equation (21). The proof is complete. 



2.10 Proof of Corollary 2 

Passing to relative error estimates requires scaling e by 1/Vk and \fk respectively. Recall 
k > 1. 

— If we invoke Theorem 1 for Go, with an absolute error e/s/k in place of e, then the number 
of terms |Ao| becomes 0((/ce~ 2 ) 1+1 /°°). Hence if we can show that \j\fk < C\\K°\\2, 
then (15) follows. 

— If we invoke Theorem 1 for G±, with an absolute error e\fk in place of e, then the 
number of terms |Ai| becomes 0(k 1+1 ^ 00 e~ 1 ^ oc + (ke~ 2 ) 1 ^ 3+1 ^°°). Hence if we can show 
that Vk < CH^sT 1 1| 2 , then (16) follows. 

— The combination of the above two error bounds would show (17), provided we can show 
that Vk > C'WK 1 — ir]K°\\2 when 77 x k. 

Therefore, it suffices to establish the lower bounds on ||if°||2, H-K" 1 H2, and H-K^ ' 1 )^. By the 
tight frame property of wave atoms, the claim for K° is exactly Jj ^ 2 \Go(s, t)\ 2 dsdt > C/k. 

Set 4>(s,t) = ||x(s) -x(t)||. As k(f>(s,t) > c , Lemma 3 implies that \G Q (s,t)\ < C^y 1 / 2 . If 
k is sufficiently large, we can restrict the integration domain to the nonempty set kcf)(s, t) > cq, 
and directly conclude. If k is not large enough for this step, the integral is still a uniformly 
continuous and positive function of k, hence uniformly bounded away from zero. 

The claim for K 1 is J [01]2 |Gi(s, t)\ 2 dsdt > k. The non-Hankel factors in the expression of 
G\ play a minor role in evaluating this lower bound; we can consider them bounded away 
from zero on a large set S± to which the integral is restricted. Lemma 3 then implies that 
for k(p(s,t) > ci and (s,t) € Si, we have \Gi(s,t)\ > k(k(f>)^ 1 / 2 . By the same reasoning as 
previously, this leads to the lower bound. 

The claim for is Jj ^ 2 |Gi(s, t) — ir/Go(s, t)\ 2 dsdt > k, with rj x k. The reasoning is here 

a little more complicated since Gi and kGo are on the same order of magnitude. The presence 
of —i, however, prevents major cancellations — and is in fact chosen for that very reason. The 
asymptotic decay of Gi — iijGq for large <f>(s, t) can be studied from the integral formulation of 
the Hankel function used throughout the Appendix for proving the three Lemmas in section 
2.1. Without entering into details, we remark that the integral factor in (56) is for z large 
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very near real-valued, with positive real part. The exponential factor e ~ mn / 2 shows that Hi 
is then almost aligned with —iHq. The particular combination G\ — ir/Go with rj > respects 
this quadrature property of Hankel functions, and produces no cancellation at all in the limit 
z — > oo. So for k<j> large enough, the claim follows; and if k<f> is not large enough, we fall back 
on an argument of uniform continuity as previously. 

3 Numerical experiments 

In this section, we provide several numerical examples to support the sparsity results of the previous 
section. The three geometric objects used in this section are displayed in Figure 1. For each object, 
the boundary curve is represented using a small number of Fourier coefficients. The last two 
examples have non-convex shapes that typically result in multiple scattering effects. We report the 
numerical results for the single layer kernel Go(s,t) in Section 3.1 and the results for the double 
layer kernel G±(s, t) in Section 3.2. We omit the results of the combined kernel Gi(s, t) — ir]Go(s, t) 
as they are almost the same as the single layer case. 




-1 -0.5 0.5 1 -1 -0.5 0.5 1 -1 -0.5 0.5 1 



(a) (b) (c) 

Figure 1: The geometric objects used in the test examples, (a): an ellipse, (b): a kite-shaped 
object, (c): a star-shaped object. 

3.1 Single layer potential 

We first study the single layer potential 

k-G (s,t) = k- l -H^\kMs)-^t)\\)\\mi 

Notice that we use k ■ Go(s,t) instead of Go(s, t), because the coupling constant 77 in the integral 
equation (1) is of order k. Therefore, k ■ Gq(s, t) is more informative when we report the value at 
which coefficients are thresholded, and the number of nonnegligible coefficients. 

For each fixed k, we construct the discrete version of the operator k ■ Go(s,t) by sampling the 
boundary curve with N = 8k quadrature points; this corresponds to about 8 points per wave- 
length in these examples. Next, we scale the values at these quadrature points with the high-order 
corrected trapezoidal quadrature rule from [20] in order to integrate the logarithmic singularity 
accurately. This quadrature rule has the appealing feature of changing the weights only locally 
close to the singularity. We then apply the two dimensional wave atom transform to compute the 
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coefficients := (kGQ,cp^). For a fixed accuracy e, we obtain the sparsest approximant K® that 
satisfies 

\\K -k \U 2M <e\\K°\\ Hfl) 

by choosing the largest possible threshold value 5 and setting the coefficients less than 5 in modulus 
to zero. Equation (15) predicts that, as a function of k, the number of wave atom coefficients 
defining K° should grow like k 1+1 /°°. 
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e = 10 -1 


e = 10" 1 - 5 


e= 10" 2 


k = 32 
k = 64 
k = 128 
k = 256 
k = 512 
= 1024 


11 / 3.26e-2 / 9.53e-2 
9 / 3.37e-2 / 1.07e-l 

10 / 3.33e-2 / 9.80e-2 
9 / 3.63e-2 / 9.79e-2 

11 / 3.61e-2 / 9.79e-2 
9 / 4.00e-2 / 9.78e-2 


18 / 9.09e-3 / 2.88e-2 

16 / 9.20e-3 / 3.39e-2 
20 / 8.46e-3 / 3.29e-2 
18 / 9.33e-3 / 2.99e-2 
20 / 9.39e-3 / 3.07e-2 

17 / 1.01e-2 / 3.12e-2 


27 / 2.60e-3 / 7.93e-3 

28 / 2.39e-3 / 1.01e-2 

32 / 2.45e-3 / 1.04e-2 
30 / 2.46e-3 / 9.66e-3 

33 / 2.55e-3 / 9.86e-3 

29 / 2.72e-3 / 9.81e-3 



Table 1: Single layer potential for the ellipse. Top left: the real part of the operator for k = 128. 
Top right: the sparsity pattern of the operator under the wave atom basis for k = 128 and e = 10~ 2 . 
Each black pixel stands for a nonnegligible coefficient. Bottom: For different combinations of k and 
e, the number of nonnegligible entries per row | Ao|/-/V, the threshold value 5, and the estimated I? 
operator norm e^2 (from left to right). 

For each example in Figure 1, we perform the test for different combinations of (k,e) with 
k = 32, 64, . . . , 1024 and e = 10 , 10 -1,5 , and 10 -2 . The numerical results for the three examples 
are summarized in Tables 1, 2, and 3, respectively. In each table, 

• The top left plot is the real part of the single layer potential in the case of k = 128. This plot 
displays coherent oscillatory patterns for which the wave atom frame is well suited. 

• The top right plot is the sparsity pattern of the operator under the wave atom basis for k = 128 
and e = 10~ 2 . Each black pixel stands for a nonnegligible coefficient. The coefficients are 
organized in a way similar to the usual ordering of 2D wave atom coefficients: each block 
contains the wave atom coefficients of a fixed frequency index (J, m), and the blocks are 
ordered such that the lowest frequency is located at the top left corner while the highest 
frequency at the bottom right corner. Within a block, the wave atom coefficients of frequency 
index (j, m) are ordered according to their spatial locations. The multiscale nature of the 
wave atom frame can be clearly seen from this plot. 
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e = 10" 1 


e = lO" 15 


e = 10- 2 


k = 32 
fc = 64 
fc = 128 
k = 256 
= 512 
k = 1024 


14 / 2.92e-2 / 9.50e-2 

15 / 2.70e-2 / 9.77e-2 
17 / 2.68e-2 / 1.04e-l 
17 / 2.70e-2 / 1.06e-l 
17 / 2.85e-2 / l.lle-1 
17 / 2.89e-2 / 1.03e-l 


25 / 8.35e-3 / 3.07e-2 
30 / 7.66e-3 / 3.45e-2 

34 / 7.49e-3 / 3.35e-2 

35 / 7.27e-3 / 3.36e-2 
35 / 7.58e-3 / 3.43e-2 
35 / 7.75e-3 / 3.33e-2 


37 / 2.69e-3 / 8.79e-3 
46 / 2.41e-3 / 9.98e-3 
53 / 2.26e-3 / 1.04e-2 
58 / 2.07e-3 / 1.05e-2 
58 / 2.09e-3 / 1.07e-2 
60 / 2.07e-3 / 1.07e-2 



Table 2: Single layer potential for the kite-shaped object. Top left: the real part of the operator 
for k = 128. Top right: the sparsity pattern of the operator under the wave atom basis for k = 128 
and e = 10 -2 . Bottom: For different combinations of k and e, \A.q\/N, 5, and e^2. 
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e = 10" 1 


e = lO" 1 - 5 


e = 10~ 2 


k = 32 
k = 64 
k = 128 
k = 256 
k = 512 
k = 1024 


14 / 2.76e-2 / 1.00e-l 
14 / 2.55e-2 / 8.97e-2 
18 / 2.39e-2 / 8.36e-2 

18 / 2.32e-2 / 9.48e-2 

19 / 2.34e-2 / 9.30e-2 
18 / 2.42e-2 / 9.50e-2 


23 / 8.74e-3 / 3.14e-2 
25 / 7.68e-3 / 2.83e-2 
33 / 7.01e-3 / 2.64e-2 
35 / 6.66e-3 / 2.90e-2 
40 / 6.19e-3 / 2.94e-2 
38 / 6.41e-3 / 3.01e-2 


32 / 2.50e-3 / 9.72e-3 
38 / 2.39e-3 / 8.28e-3 
49 / 2.19e-3 / 8.06e-3 
56 / 1.98e-3 / 9.33e-3 
66 / 1.84e-3 / 9.30e-3 
66 / 1.77e-3 / 9.41e-3 



Table 3: Single layer potential for the star-shaped object. Top left: the real part of the operator 
for k = 128. Top right: the sparsity pattern of the operator under the wave atom basis for k = 128 
and e = 10 -2 . Bottom: For different combinations of k and e, \Aq\/N, 5, and e^2. 
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• The table at the bottom gives, for different combinations of k and e, the number of non- 
negligible coefficients per row |Ao|/iV , the threshold value 5 (coefficients below this value in 
modulus are put to zero) and the L 2 -to-L 2 norm operator error e L 2 estimated using random 
test functions. 

In these tables, the number of significant coefficients per row Aq/N grows very slowly as k 
doubles and reaches a constant level for large values of k. This match well with the theoretical 
analysis in Section 2. The threshold value 5 remains roughly at a constant level as k grows, which 
is quite different from the results obtained using wavelet packet bases [13, 14, 17, 18] where the 
threshold value in general decreases as k grows. The estimated L 2 -to-L 2 operator error e^2 is very 
close to the prescribed accuracy e in all cases. This indicates that, in order to get an approximation 
within accuracy e in operator norm, one can simply truncate the nonstandard form of the operator 
in the wave atom frame with the same accuracy. 



3.2 Double layer potential 

We now consider the double layer potential 

Gi(a,t) = jH^ikMs) -x(t)||) ii^I^n -"xM ll*(*)ll- 

For each fixed k, the discrete version of G\(s,t) is constructed by sampling at N = 8k points 
and using trapezoidal quadrature rule. The coefficients K\ := (G\,ip^) are calculated using the 
two dimensional wave atom transform and the approximant is constructed in the same way as 
the single layer potential case. 

The results of the double layer potentials for the three examples are summarized in Tables 4, 5, 
and 6, respectively. These results are qualitatively similar to the ones of the singular layer potential. 
However, the coefficients of the double layer potential exhibit better sparsity pattern for the simple 
reason that the double layer potential operator has a singularity much weaker than logarithmic 
along the diagonal (where s = t) for objects with smooth boundary. Therefore, for a fixed accuracy 
e, the number of wave atoms required along the diagonal for the double layer potential is smaller 
than the number for the singular layer potential. This is clearly shown in the sparsity pattern plots 
in Tables 4, 5, and 6. 



A Additional proofs 

Proof of Lemma 1. 

Following Watson's treatise [26], the Hankel function can be expressed by complex contour 
integration as 

where —tt/2 < [3 < tt/2. For us, z is real and positive, and we take (3 = for simplicity. 

Let us first treat the case m = (no differentiations) and n > 0. We can use the simple bound 

1 + — 2 < C n 1+ - 3) 

2z \ z J 

to see that the integral, in absolute value, is majorized by C n (1 + z~ n+l l 2 ). Hence the Hankel 
function itself is bounded by C n (z -1 / 2 + z~ n ). This establishes the first two expressions in (22) in 
the case m = 0. 
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e = lO" 1 


e = 10~ L5 


e = 10- 2 


k = 32 
k = 64 
ft; = 128 
ft; = 256 
k = 512 
fc = 1024 


10 / 1.39e-2 / 1.29e-l 
9 / 1.40e-2 / 1.00e-l 

11 / 1.21e-2 / 9.50e-2 

10 / 1.28e-2 / 9.80e-2 
13 / 1.16e-2 / 1.00e-l 

11 / 1.21e-2 / 9.98e-2 


16 / 4.61e-3 / 3.65e-2 
15 / 4.41e-3 / 3.23e-2 
20 / 3.74e-3 / 3.00e-2 
18 / 3.70e-3 / 3.28e-2 
22 / 3.42e-3 / 3.24e-2 
20 / 3.52e-3 / 3.29e-2 


21 / 1.57e-3 / 1.196-2 
21 / 1.34e-3 / 1.01e-2 
28 / 1.12e-3 / 9.94e-3 
28 / 1.07e-3 / 1.02e-2 
33 / 9.77e-4 / 1.00e-2 
30 / 9.83e-4 / 1.01e-2 



Table 4: Double layer potential for the ellipse. Top left: the real part of the operator for k = 128. 
Top right: the sparsity pattern of the operator under the wave atom basis for k = 128 and e = 10 -2 . 
Bottom: For different combinations of k and e, |Ai|/iV, 5, and e L 2. 
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e = 10 _1 


e = 10 -1 - 6 


e= 10" 2 


k = 32 
k = 64 
k = 128 
k = 256 
k = 512 
k = 1024 


12 / 1.64e-2 / 9.16e-2 

13 / 1.40e-2 / 1.02e-l 
16 / 1.27e-2 / 9.11e-2 
16 / 1.22e-2 / 9.56e-2 
18 / 1.20e-2 / 9.98e-2 
18 / 1.17e-2 / 9.57e-2 


19 / 5.16e-3 / 2.82e-2 
23 / 4.37e-3 / 3.18e-2 
29 / 3.82e-3 / 2.80e-2 
31 / 3.46e-3 / 3.23e-2 
33 / 3.42e-3 / 3.09e-2 
35 / 3.26e-3 / 3.02e-2 


27 / 1.64e-3 / 8.71e-3 
34 / 1.31e-3 / 1.10e-2 
43 / 1.18e-3 / 9.27e-3 
49 / 1.02e-3 / 1.01e-2 
52 / 9.71e-4 / 9.81e-3 
55 / 9.07e-4 / 1.03e-2 



Table 5: Double layer potential for the kite-shaped object. Top left: the real part of the operator 
for k = 128. Top right: the sparsity pattern of the operator under the wave atom basis for k = 128 
and e = 10 -2 . Bottom: For different combinations of k and e, |Ai|/iV, 6, and e^2. 
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e = lO" 1 


e = lO" 15 


e = lO' 2 


fc = 32 
fc = 64 
fe = 128 
= 256 
k = 512 
k = 1024 


13 / 1.43e-2 / 1.20e-l 

14 / 1.25e-2 / 9.20e-2 

19 / 1.08e-2 / 1.08e-l 

20 / 1.00e-2 / 1.02e-l 
23 / 8.99e-3 / 1.01e-l 
22 / 9.06e-3 / 1.04e-l 


20 / 5.15e-3 / 3.31e-2 
23 / 4.23e-3 / 2.80e-2 
32 / 3.42e-3 / 3.27e-2 
36 / 3.08e-3 / 3.13e-2 
44 / 2.66e-3 / 3.24e-2 
43 / 2.54e-3 / 3.21e-2 


25 / 1.62e-3 / 1.15e-2 
31 / 1.43e-3 / 9.03e-3 
44 / 1.19e-3 / 1.10e-2 
52 / 9.77e-4 / 1.07e-2 
66 / 8.22e-4 / 1.01e-2 
69 / 7.46e-4 / 9.46e-3 



Table 6: Double layer potential for the star-shaped object. Top left: the real part of the operator 
for k = 128. Top right: the sparsity pattern of the operator under the wave atom basis for k = 128 
and e = 10 -2 . Bottom: For different combinations of k and e, |Ai|/iV, 6, and e^2. 



The case m = n = is treated a little differently because the integrand in (56) develops a 1/u 
singularity near the origin as z — > 0. We have 

hence the integral in (56) is bounded in modulus by a constant times 

pz roc 

/ e- u u- l / 2 du + z 1 / 2 e- u u^du < C {z 1 / 2 + z l / 2 \\ogz\). 

JO Jz 

With the z -1 / 2 factor from (56), the resulting bound is C (1 + | logz|) as desired (third equation). 
When z > 1, we can improve this to C(l + z 1 / 2 e- z ) < C, which gives the first equation when 
m = n = 0. 

For the case m > 0, it suffices to apply Leibniz's rule inductively and observe that each derivative 
produces a factor x~ l without changing the power of k. In particular, 

• With ti^O, 

_ [(kx)~ a ] = -a{kx)- a -, 
hence the power of k is preserved and one negative power of x is created. 

• Derivatives acting on isolated negative powers of x also produce ai" 1 factor without affecting 
the dependence on k. 
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• As for the x-dependence under the integral sign, with a / 0, we arrange the factors as 



d I ( iu V 
dx~ [ V + 2fcx~y 



Q: 



in / iu \ 
2kx V + 2fccy) 



1 + ^ U \ a ^ 
2kx J x 



The factor in square brackets is bounded by 1 in modulus, hence the dependence on k is 
unchanged. The factor 1/x is the only modification in the dependence on x. Subsequent 
differentiations will only act on factors that we have already treated above: powers of kx, 
powers of x, and powers of 1 + This finishes the proof. 



Proof of Lemma 2. 

As previously, we use the integral formulation to get 



f'OO 

xH[ 1) (x) = f(x) J e-V/ 2 (x + i^j 



du, 



where /(x) is the exponential factor, and already obeys \f^ n \x)\ < C for all n > 0. We denote the 
integral factor by I(x); its derivatives are 



jW(x) = C n J™ e^u 1 ' 2 (x + ^ n du, 



where C n is a numerical constant. In a manner analogous to the proof of Lemma 1, we can bound 



\x + i||2" n = ( x 2 + ^j- ) < C n (max(x,n))2 n 



^ \ 4 2 



It follows that 



T (•'£ poo 

\I (n) (x)\<C n / e- u u 1/2 x 1/2 - n du+ e~ u u l - n du . 



In the first term we can use e~ u < 1 and bound the integral by a constant times x 2 " ra . The 
integrand of the second term has a singularity near u = that becomes more severe as n increases; 
this term is bounded by 0(1) if n = or 1, by 0(1 + | log(x)|) if n = 2, and by 0(x 2 ~ n ) if n > 2. 

Proof of Lemma 3. 

Consider equation (56) again, and take z real. For large values of z, the factor (1 + iuj2z) n ~ x l 2 
is close to 1; more precisely, it is easy to show that for each n > 0, there exists c n > 0, d n > for 
which 

zu \ u 

1 + — -l\<d n -, if u<c n z. 

2z J z 

We can insert this estimate in (56) and split the integral into two parts to obtain 
|^)(.)|r(n-l/2)(^) 1/2 - \£ Z e- u u-^ 2 du\> 



rzc„ /-oo 

d n e ^ u ~ l / 2 -du-C e- u u x l 2 du. 

JO Z Jzc n 



The constant C in the last term comes from equation (57). The first term in the right-hand side is a 
0(z~ r ), and the second term is a 0{e~ z ). At the expense of possibly choosing increasing the value 
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of c n , the second term in the left-hand side can manifestly be made to dominate the contribution 
of the right hand side, proving the lemma. 

Proof of Lemma 6. 

We start by writing 

(x(t))^-r = (x(0) ± - X(t l + f (t) "^ ) . 

ll x 0) - x (s)ll 

By Taylor's theorem, |x(t) + crx(i) — x(s)| < Ca 2 , hence Kx(t))- 1 - • r| < Co. Derivatives are then 
treated by induction; recall that ^ = ^ + and ^ = ^ — ^; 

• any number of r or cr derivatives acting on (x(t))- 1 leave it a O(l); 

• r derivatives acting on x(t) + crx(t) — x(s) leave it a 0(<r 2 ) while each u derivative removes 
an order of a; 

• r derivatives acting on ||x(t)— x(s) || _m leave it a 0(cr~ m ) (m is generic) while each a derivative 
removes an order of a. 

This shows (50) and (51). 
Proof of Lemma 7. 

Some cancellations will need to be quantified in this proof, that were not a concern in the 
justification of previous coarser estimates like Lemma 6. 

Without loss of generality, assume that x(i) = (0, 0), n^u) = (0, 1), and that we have performed 
a change of variables such that the curve is parametrized as the graph x(s) = (s,/(s)) of some 
function / € C°° obeying < Cs 2 . This latter change of variables would contribute a bounded 

multiplicative factor that would not compromise the overall estimate. 

By symmetry, if (54) is true for d/ds derivatives, then it will be true for d/dt derivatives as 
well. Without loss of generality let s > 0. Then we have 



0(5, t) = Ms) - x(t) || = Vs 2 + f 2 (s) =syjl + ^-. 

Since \f(s)\ < Cs 2 and C°°, the ratio f 2 (s)/s 2 is also bounded for s < 1, and of class C°°. Being 
a composition of C°° functions, the whole factor + f 2 /s 2 is therefore also of class C°°, which 
proves (54). 

As for (55) with d/ds derivatives, we can write (x(t))- 1 = (0, 1) and 
{±{t)) ± x(t)-x(s) _ -f(s) _ ff(s) 



|x(t)-x(s)||2 S 2 + f 2 ( S ) - V * 2 J 1 + 

(The l/||x(t)|| does not pose a problem since it is C°°.) Again, since \f(s)\ < Cs 2 , both ratios 
f(s)/s 2 and f 2 (s)/s 2 are themselves bounded and of class C°°. The factor — is the composition 

s z 

of two C°° functions, hence also of class C°°. 

The symmetry argument is not entirely straightforward for justifying (55) with d/dt derivatives. 
Symmetry s <-> t only allows to conclude that 



|x(t)-x( S )|P 



< Cm) 
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where (x(s)) -1 " appears in place of the desired (x^))- 1 . Hence it suffices to show that the d/dt, or 
equivalently the d/ds derivatives of 



(x(t))- 1 - (X( S ))" 



x(t) - x(s) 



||x(t)-x( S )P 

stay bounded. Using our frame in which the curve is a graph, we find 

(x(t))^ - (i(s)) ± = (0,1) - (-/'(,), 1) = (/'(*), 0) 

As a result, 



(kit)) 1 - - (x( S ))- 



x(t)-x(s) -g/^g) 

|x(t)-x(s)||2 S 2 + /2(s) 



1 + /2{S) 
1 H 7T~ 



Since |/'(s)| < s, we are again in presence of a combination of C°° functions that stays infinitely 
differentiable. 
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